import scanpy as sc, anndata as ad, numpy as np, pandas as pd
from scipy import sparse
from anndata import AnnData
import warnings
import socket
import plotly.express as px
from matplotlib import pylab
import sys
import yaml
import os
import matplotlib.pyplot as plt
import scvelo as scv
from pybiomart import Server
import rpy2
from pandas.api.types import CategoricalDtype
import rpy2.robjects as ro
import seaborn as sns
from matplotlib_venn import venn2, venn2_circles, venn2_unweighted
from matplotlib_venn import venn3, venn3_circles
from matplotlib import pyplot as plt
import plotly.graph_objects as go
import ipyparams
warnings.filterwarnings('ignore')
/group/testa/Users/davide.castaldi/Polaroid_publication_revisions/.local/lib/python3.8/site-packages/setuptools_scm/_integration/setuptools.py:30: RuntimeWarning: ERROR: setuptools==56.1.0 is used in combination with setuptools_scm>=8.x Your build configuration is incomplete and previously worked by accident! setuptools_scm requires setuptools>=61 Suggested workaround if applicable: - migrating from the deprecated setup_requires mechanism to pep517/518 and using a pyproject.toml to declare build dependencies which are reliably pre-installed before running the build tools warnings.warn(
outdir="./outdir"
FinaLeaf="/Neurons"
output_basename = "09A_Neurons_pBulk_bySegment.DEA"
mappingDict = {}
categoriesOrder = ["piece1","piece2","piece3","distal","medial","proximal"]
groupingCovariate = "group"
PseudooReplicates_per_group = 10
totalPath = outdir+FinaLeaf+"/5A_Neurons_pBulk.bySegment."+str(PseudooReplicates_per_group)+"PRs.tsv"
adataPath = outdir+FinaLeaf+"/5A_Neurons_pBulk.bySegment."+str(PseudooReplicates_per_group)+"PRs.h5ad"
if not os.path.exists('./figures'):
os.makedirs('./figures')
if not os.path.exists('./tables'):
os.makedirs('./tables')
%load_ext rpy2.ipython
%%R
library(edgeR)
library(org.Hs.eg.db)
library(AnnotationDbi)
library(stats)
library(topGO)
R[write to console]: Loading required package: limma
R[write to console]: Loading required package: AnnotationDbi
R[write to console]: Loading required package: stats4
R[write to console]: Loading required package: BiocGenerics
R[write to console]: Loading required package: parallel
R[write to console]:
Attaching package: ‘BiocGenerics’
R[write to console]: The following objects are masked from ‘package:parallel’:
clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
clusterExport, clusterMap, parApply, parCapply, parLapply,
parLapplyLB, parRapply, parSapply, parSapplyLB
R[write to console]: The following object is masked from ‘package:limma’:
plotMA
R[write to console]: The following objects are masked from ‘package:stats’:
IQR, mad, sd, var, xtabs
R[write to console]: The following objects are masked from ‘package:base’:
anyDuplicated, append, as.data.frame, basename, cbind, colnames,
dirname, do.call, duplicated, eval, evalq, Filter, Find, get, grep,
grepl, intersect, is.unsorted, lapply, Map, mapply, match, mget,
order, paste, pmax, pmax.int, pmin, pmin.int, Position, rank,
rbind, Reduce, rownames, sapply, setdiff, sort, table, tapply,
union, unique, unsplit, which.max, which.min
R[write to console]: Loading required package: Biobase
R[write to console]: Welcome to Bioconductor
Vignettes contain introductory material; view with
'browseVignettes()'. To cite Bioconductor, see
'citation("Biobase")', and for packages 'citation("pkgname")'.
R[write to console]: Loading required package: IRanges
R[write to console]: Loading required package: S4Vectors
R[write to console]:
Attaching package: ‘S4Vectors’
R[write to console]: The following object is masked from ‘package:base’:
expand.grid
R[write to console]:
R[write to console]: Loading required package: graph
R[write to console]: Loading required package: GO.db
R[write to console]:
R[write to console]: Loading required package: SparseM
R[write to console]:
Attaching package: ‘SparseM’
R[write to console]: The following object is masked from ‘package:base’:
backsolve
R[write to console]:
groupGOTerms: GOBPTerm, GOMFTerm, GOCCTerm environments built.
R[write to console]:
Attaching package: ‘topGO’
R[write to console]: The following object is masked from ‘package:IRanges’:
members
%matplotlib inline
sc.settings.verbosity = 3 # verbosity: errors (0), warnings (1), info (2), hints (3)
sc.logging.print_header()
sc.settings.set_figure_params(dpi=50, facecolor='white')
pylab.rcParams['figure.figsize'] = (10, 10)
scanpy==1.8.0 anndata==0.8.0 umap==0.4.6 numpy==1.22.2 scipy==1.10.1 pandas==1.2.3 scikit-learn==0.24.2 statsmodels==0.13.1 python-igraph==0.9.8 louvain==0.7.1 pynndescent==0.5.5
adata = sc.read_h5ad(adataPath)
adata.obs[groupingCovariate] = adata.obs[groupingCovariate].astype(CategoricalDtype(categories=categoriesOrder, ordered=True))
total = pd.read_csv(totalPath, sep="\t", index_col=0)
#adata = adata[adata.obs["group"].isin(RelevantAreas)]
total = total[adata.obs_names]
edgeR_topTags = ro.r['topTags']
edgeR_glmLRT = ro.r['glmLRT']
edgeR_glmQLFTest = ro.r['glmQLFTest']
as_data_frame = ro.r['as.data.frame']
rprint = rpy2.robjects.globalenv.find("print")
DEGs = {}
TestCov = groupingCovariate
adata.obs.group.value_counts()
piece1 10 piece2 10 piece3 10 distal 10 medial 10 proximal 10 Name: group, dtype: int64
obs = adata.obs[[TestCov,"pseudoreplicate"]].copy()
obs["pseudoreplicate"] = obs.index.tolist()
obs["TestCov"] = obs[TestCov].astype(str)
#obs["TestCov"] = obs[TestCov]
totalRelevant = total[obs.index]
totalRelevant.head()
| distal_0 | distal_1 | distal_2 | distal_3 | distal_4 | distal_5 | distal_6 | distal_7 | distal_8 | distal_9 | ... | proximal_0 | proximal_1 | proximal_2 | proximal_3 | proximal_4 | proximal_5 | proximal_6 | proximal_7 | proximal_8 | proximal_9 | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| MIR1302-2HG | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.0 | 0.0 | ... | 0.564133 | 0.000000 | 0.000000 | 0.00000 | 0.000000 | 0.000000 | 0.0 | 0.00000 | 0.000000 | 0.000000 |
| FAM138A | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.0 | 0.0 | ... | 0.000000 | 0.000000 | 0.000000 | 0.00000 | 0.000000 | 0.000000 | 0.0 | 0.00000 | 0.000000 | 0.000000 |
| OR4F5 | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.0 | 0.0 | ... | 0.000000 | 0.000000 | 0.000000 | 0.00000 | 0.000000 | 0.000000 | 0.0 | 0.00000 | 0.000000 | 0.000000 |
| AL627309.1 | 0.0 | 0.0 | 0.0 | 1.823126 | 0.924461 | 0.903551 | 2.596238 | 0.845117 | 0.0 | 0.0 | ... | 1.213359 | 2.400798 | 1.462998 | 2.86158 | 1.283984 | 0.706691 | 0.0 | 4.56084 | 3.235815 | 0.613853 |
| AL627309.3 | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.0 | 0.0 | ... | 0.000000 | 0.000000 | 0.000000 | 0.00000 | 0.000000 | 0.000000 | 0.0 | 0.00000 | 0.000000 | 0.000000 |
5 rows × 60 columns
minCounts = totalRelevant[totalRelevant > 0].min(axis = 1).min()
#Universe 1 count in at least 5% of samples per group
groupThreshold = round(pd.get_dummies(obs[TestCov]).T.sum(axis = 1) *0.05)
bolMatrix = (np.matrix(np.dot(pd.get_dummies(obs[TestCov]).T, (totalRelevant > 0).astype(int).T)) - np.matrix(groupThreshold).T > 0)
bolVect = (bolMatrix.sum(axis = 0) >= 1).A1
totalRelevant = totalRelevant.loc[bolVect]
len(totalRelevant)
27149
universe = totalRelevant.index.tolist()
levelsToMap = adata.obs[groupingCovariate].cat.categories.tolist()
levelsToMap
['piece1', 'piece2', 'piece3', 'distal', 'medial', 'proximal']
%%R -i obs -i totalRelevant -i levelsToMap -o dds
mmVector<- factor(obs$TestCov, levels = c(levelsToMap))
mm <- model.matrix(~mmVector )
row.names(mm) <- colnames(totalRelevant)
dds <- DGEList(totalRelevant, group = obs$TestCov, genes = rownames(totalRelevant))
#Ennotate with entrez genes
anno <- select(org.Hs.eg.db, keys=rownames(dds$counts),columns=c("ENTREZID","SYMBOL"),keytype="SYMBOL")
colnames(anno) <- c("genes","entrez")
anno <- anno[!duplicated(anno$genes, fromLast=T), ]
dds$genes$entrez <- merge(x = data.frame(dds$genes), y = anno, by = "genes", all.x = TRUE)$entrez
dim(dds)
R[write to console]: 'select()' returned 1:many mapping between keys and columns
[1] 27149 60
%%R
mm
(Intercept) mmVectorpiece2 mmVectorpiece3 mmVectordistal
distal_0 1 0 0 1
distal_1 1 0 0 1
distal_2 1 0 0 1
distal_3 1 0 0 1
distal_4 1 0 0 1
distal_5 1 0 0 1
distal_6 1 0 0 1
distal_7 1 0 0 1
distal_8 1 0 0 1
distal_9 1 0 0 1
medial_0 1 0 0 0
medial_1 1 0 0 0
medial_2 1 0 0 0
medial_3 1 0 0 0
medial_4 1 0 0 0
medial_5 1 0 0 0
medial_6 1 0 0 0
medial_7 1 0 0 0
medial_8 1 0 0 0
medial_9 1 0 0 0
piece1_0 1 0 0 0
piece1_1 1 0 0 0
piece1_2 1 0 0 0
piece1_3 1 0 0 0
piece1_4 1 0 0 0
piece1_5 1 0 0 0
piece1_6 1 0 0 0
piece1_7 1 0 0 0
piece1_8 1 0 0 0
piece1_9 1 0 0 0
piece2_0 1 1 0 0
piece2_1 1 1 0 0
piece2_2 1 1 0 0
piece2_3 1 1 0 0
piece2_4 1 1 0 0
piece2_5 1 1 0 0
piece2_6 1 1 0 0
piece2_7 1 1 0 0
piece2_8 1 1 0 0
piece2_9 1 1 0 0
piece3_0 1 0 1 0
piece3_1 1 0 1 0
piece3_2 1 0 1 0
piece3_3 1 0 1 0
piece3_4 1 0 1 0
piece3_5 1 0 1 0
piece3_6 1 0 1 0
piece3_7 1 0 1 0
piece3_8 1 0 1 0
piece3_9 1 0 1 0
proximal_0 1 0 0 0
proximal_1 1 0 0 0
proximal_2 1 0 0 0
proximal_3 1 0 0 0
proximal_4 1 0 0 0
proximal_5 1 0 0 0
proximal_6 1 0 0 0
proximal_7 1 0 0 0
proximal_8 1 0 0 0
proximal_9 1 0 0 0
mmVectormedial mmVectorproximal
distal_0 0 0
distal_1 0 0
distal_2 0 0
distal_3 0 0
distal_4 0 0
distal_5 0 0
distal_6 0 0
distal_7 0 0
distal_8 0 0
distal_9 0 0
medial_0 1 0
medial_1 1 0
medial_2 1 0
medial_3 1 0
medial_4 1 0
medial_5 1 0
medial_6 1 0
medial_7 1 0
medial_8 1 0
medial_9 1 0
piece1_0 0 0
piece1_1 0 0
piece1_2 0 0
piece1_3 0 0
piece1_4 0 0
piece1_5 0 0
piece1_6 0 0
piece1_7 0 0
piece1_8 0 0
piece1_9 0 0
piece2_0 0 0
piece2_1 0 0
piece2_2 0 0
piece2_3 0 0
piece2_4 0 0
piece2_5 0 0
piece2_6 0 0
piece2_7 0 0
piece2_8 0 0
piece2_9 0 0
piece3_0 0 0
piece3_1 0 0
piece3_2 0 0
piece3_3 0 0
piece3_4 0 0
piece3_5 0 0
piece3_6 0 0
piece3_7 0 0
piece3_8 0 0
piece3_9 0 0
proximal_0 0 1
proximal_1 0 1
proximal_2 0 1
proximal_3 0 1
proximal_4 0 1
proximal_5 0 1
proximal_6 0 1
proximal_7 0 1
proximal_8 0 1
proximal_9 0 1
attr(,"assign")
[1] 0 1 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$mmVector
[1] "contr.treatment"
%%R -i minCounts
#keep <- filterByExpr(dds, min.count = minCounts, min.total.count = minCounts*5)
#dds <- dds[keep,,keep.lib.sizes=FALSE]
dim(dds)
[1] 27149 60
%%R
dds <- calcNormFactors(dds)
dds <- estimateGLMRobustDisp(dds, mm)
fit <- glmFit(dds, mm)
levelsToMap
['piece1', 'piece2', 'piece3', 'distal', 'medial', 'proximal']
%%R -o piece2_vs_piece1 -o piece3_vs_piece1 -o distal_vs_piece1 -o medial_vs_piece1 -o proximal_vs_piece1 -o proximal_vs_distal -o medial_vs_distal -o proximal_vs_medial
#-o Mid_vs_early -o Late_vs_mid
topn=25000
piece2_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=2), adjust.method = "fdr", n = topn, sort.by = "logFC"))
piece3_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=3), adjust.method = "fdr", n = topn, sort.by = "logFC"))
distal_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=4), adjust.method = "fdr", n = topn, sort.by = "logFC"))
medial_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=5), adjust.method = "fdr", n = topn, sort.by = "logFC"))
proximal_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=6), adjust.method = "fdr", n = topn, sort.by = "logFC"))
proximal_vs_distal = as.data.frame(topTags(glmLRT(fit, contrast=c(0,0,0,-1,0,1)), adjust.method = "fdr", n = topn, sort.by = "logFC"))
medial_vs_distal = as.data.frame(topTags(glmLRT(fit, contrast=c(0,0,0,-1,1,0)), adjust.method = "fdr", n = topn, sort.by = "logFC"))
proximal_vs_medial = as.data.frame(topTags(glmLRT(fit, contrast=c(0,0,0,0,-1,1)), adjust.method = "fdr", n = topn, sort.by = "logFC"))
DEGsDict = {}
DEGsDict = {"piece2_vs_piece1":piece2_vs_piece1,
"piece3_vs_piece1":piece3_vs_piece1,
"distal_vs_piece1":distal_vs_piece1,
"medial_vs_piece1":medial_vs_piece1,
"proximal_vs_piece1":proximal_vs_piece1,
"proximal_vs_distal":proximal_vs_distal,
"medial_vs_distal":medial_vs_distal,
"proximal_vs_medial":proximal_vs_medial,
}
LOGCFTHRESHOLD = 1.5
QTHRESH = 0.05
for k in list(DEGsDict.keys()):
DEGsDict[k]["-logPVal"] = -np.log10(DEGsDict[k].PValue)
DEGsDict[k]["significant"] = "notSignificant"
DEGsDict[k].loc[(DEGsDict[k]["logFC"] < -LOGCFTHRESHOLD) & (DEGsDict[k]["FDR"] < QTHRESH),"significant"] = "Down"
DEGsDict[k].loc[(DEGsDict[k]["logFC"] > LOGCFTHRESHOLD) & (DEGsDict[k]["FDR"] < QTHRESH),"significant"] = "Up"
color_discrete_map = {'notSignificant': 'white', 'Up': 'red', 'Down': 'blue'}
fig = px.scatter(DEGsDict[k], x="logFC", y="-logPVal",
#hover_data=DEGsDict[k][["genes","knownMarker"]],
hover_data=DEGsDict[k][["genes"]],
width=600, height=600,
color = "significant", color_discrete_map=color_discrete_map,
#symbol="knownMarker",
title=k,
template="simple_white")
fig.update_traces(mode='markers', marker_line_width=1, marker_size=10)
fig.add_hline(y=DEGsDict[k].loc[DEGsDict[k]["FDR"] - 0.01 < 0,"-logPVal" ].min(),
line_color="black",line_width=2, line_dash="dash")
fig.add_vline(x=LOGCFTHRESHOLD,
line_color="black",line_width=2, line_dash="dash")
fig.add_vline(x=-LOGCFTHRESHOLD,
line_color="black",line_width=2, line_dash="dash")
fig.show()
for k in DEGsDict.keys():
DEGsDict[k][DEGsDict[k].significant != "notSignificant" ].sort_values("logFC", ascending = False).to_csv(outdir+FinaLeaf+"/"+output_basename+"."+ k + ".tsv", sep = "\t")
DEGsDict[k][DEGsDict[k].significant != "notSignificant" ].sort_values("logFC", ascending = False).to_excel(outdir+FinaLeaf+"/"+output_basename+"."+ k + ".xlsx")
proximal_vs_distalDF = DEGsDict["proximal_vs_distal"]
proximal_vs_distalDF["Contrast"] = "proximal_vs_distal"
medial_vs_distalDF = DEGsDict["medial_vs_distal"]
medial_vs_distalDF["Contrast"] = "medial_vs_distal"
OverallContrsts = pd.concat([proximal_vs_distalDF,medial_vs_distalDF], ignore_index=True)
#OverallContrsts = OverallContrsts[OverallContrsts["significant"] != "notSignificant"]
gp = OverallContrsts.groupby(["genes","significant"], as_index=False).size()
UpAtl1 = gp[(gp["significant"] == "Up") & (gp["size"] >=1)].genes.tolist()
DownAtl1 = gp[(gp["significant"] == "Down") & (gp["size"] >=1)].genes.tolist()
OverallContrsts = OverallContrsts[OverallContrsts.genes.isin(DownAtl1+UpAtl1)]
# OverallContrsts = pd.concat([OverallContrsts,pd.DataFrame({"genes":OverallContrsts.genes.unique().tolist(),
# "logFC":0,"significant":"notSignificant",
# "Contrast":"distal"})])
OverallContrsts["significant"] = pd.Categorical(OverallContrsts.significant, categories = ['notSignificant','Up', 'Down'],ordered=True)
OverallContrsts["Contrast"] = pd.Categorical(OverallContrsts.Contrast, categories = ['proximal_vs_distal','medial_vs_distal','distal'],ordered=True)
OverallContrsts["ABSlogFC"] = abs(OverallContrsts["logFC"])
#OverallContrsts = OverallContrsts[OverallContrsts["significant"] != "notSignificant"]
#sns.catplot(data=OverallContrsts, x="Contrast", y="logFC", hue="logFC")
fig, ax = plt.subplots(1,1, figsize=(20,15))
sns.axes_style("ticks")
sns.scatterplot(data=OverallContrsts.sort_values(['Contrast','significant']), palette=color_discrete_map,
x="Contrast", y="logFC", size="-logPVal", hue="significant", ax=ax, linewidth=.2,sizes=(100, 700),
edgecolor="black")
sns.lineplot(data=OverallContrsts, x="Contrast", y="logFC",units="genes", estimator=None, alpha=.05, color="grey", ax = ax, zorder=-1)
sns.despine()
ax.legend(bbox_to_anchor=(1.22, 0.8), prop={'size': 30})
# Plot labels
nGenes = 10
# Get pixels size
bbox = ax.get_window_extent().transformed(fig.dpi_scale_trans.inverted())
widthPXL, heightPXL = bbox.width, bbox.height
heightPXL *= fig.dpi
width = 1
height = 1
QuadrantFractionDictH = dict(zip(["Down","notSignificant","Up"],
[1,0,-1]))
offset = lambda p: transforms.ScaledTranslation(p/72.,0, plt.gcf().dpi_scale_trans)
trans = plt.gca().transData
ContrastTuple = (OverallContrsts.loc[OverallContrsts['significant'] != "notSignificant",["Contrast",'significant']].drop_duplicates())
ContrastTuple = list(zip(ContrastTuple.Contrast, ContrastTuple.significant))
for i in ContrastTuple:
PointsText = OverallContrsts[(OverallContrsts["Contrast"] == i[0]) & (OverallContrsts["significant"] == i[1])].sort_values("ABSlogFC", ascending=False).head(nGenes)
PointsText["TextRank"] = (-PointsText["ABSlogFC"]).rank()-1
PointsText["yText"] = [QuadrantFractionDictH[i] for i in PointsText["significant"]]
PointsText["yText"] = PointsText["TextRank"]*PointsText["yText"]
PointsText = PointsText.sort_values("TextRank", ascending=True)
PointsText["FCdelta"] = PointsText["logFC"] - PointsText["logFC"].values[0]
PointsText["yText2"] = PointsText["yText"] - PointsText["FCdelta"]
PointsText["xTextSwitch"] = np.where(PointsText["Contrast"] == "medial_vs_distal",-1,1)
PointsText["xAnchorSwitch"] = np.where(PointsText["Contrast"] == "medial_vs_distal","right","left")
PointsText["xAnchorSwitchArrow"] = np.where(PointsText["Contrast"] == "medial_vs_distal",1,0)
#point["xAnchorSwitch"]
for i, point in PointsText.iterrows():
ax.annotate(str(point['genes']), va="center", ha=point["xAnchorSwitch"],
xy=(point["Contrast"], point['logFC']),
xytext=(200*point["xTextSwitch"], (point["yText2"]*heightPXL/26)),
textcoords='offset points',fontsize=25,
arrowprops=dict(facecolor='black', arrowstyle='-',relpos=(point["xAnchorSwitchArrow"],0.5)))
plt.title("DEGs: {} Qval:{} MinlogFC:{}".format(output_basename.split("_")[1], QTHRESH,LOGCFTHRESHOLD), fontsize=30)
plt.tick_params(axis='both', labelsize=30)
plt.xlabel('Contrast', fontsize=40)
plt.ylabel('logFC', fontsize=40)
fig.savefig("./figures/Volcanos-{}_Qval-{}_MinlogFC-{}.png".format(output_basename.split("_")[1], QTHRESH,LOGCFTHRESHOLD), bbox_inches='tight',dpi = 400)
proximal_vs_distalDF = DEGsDict["proximal_vs_distal"]
proximal_vs_distalDF["Contrast"] = "proximal_vs_distal"
medial_vs_distalDF = DEGsDict["medial_vs_distal"]
medial_vs_distalDF["Contrast"] = "medial_vs_distal"
OverallContrsts = pd.concat([proximal_vs_distalDF,medial_vs_distalDF], ignore_index=True)
#OverallContrsts = OverallContrsts[OverallContrsts["significant"] != "notSignificant"]
# OverallContrsts = pd.concat([OverallContrsts,pd.DataFrame({"genes":OverallContrsts.genes.unique().tolist(),
# "logFC":0,"significant":"notSignificant",
# "Contrast":"distal"})])
OverallContrsts["significant"] = pd.Categorical(OverallContrsts.significant, categories = ['notSignificant','Up', 'Down'],ordered=True)
OverallContrsts["Contrast"] = pd.Categorical(OverallContrsts.Contrast, categories = ['proximal_vs_distal','medial_vs_distal','distal'],ordered=True)
OverallContrsts["ABSlogFC"] = abs(OverallContrsts["logFC"])
OverallContrstsSS = OverallContrsts[(OverallContrsts["significant"] != "notSignificant") & (OverallContrsts["Contrast"] != "distal")]
OverallContrstsSS["Contrast"] = OverallContrstsSS["Contrast"].astype("string")
OverallContrstsSS["significant"] = OverallContrstsSS["significant"].astype("string")
OverallContrstsSS = OverallContrstsSS.groupby(["Contrast","significant"], as_index=False).size()
OverallContrstsSS["Contrast"] = pd.Categorical(OverallContrstsSS.Contrast, categories = ['proximal_vs_distal','medial_vs_distal'],ordered=True)
import matplotlib.patches as mpatches
sns.axes_style("ticks",{'axes.grid' : False})
fig, (ax1, ax2) = plt.subplots(ncols=2, figsize=(20, 10), sharey=True, gridspec_kw={'wspace': 0})
sns.barplot(data=OverallContrstsSS[OverallContrstsSS['significant'] == 'Up'], x='size', y='Contrast',color="red",
ci=False, orient='horizontal', dodge=True, ax=ax2)
ax2.yaxis.set_label_position('right')
ax2.tick_params(axis='y', labelright=True, right=True)
ax2.set_title(' '+'Upregulated Genes', fontsize=40)
ax2.set_xlabel('number of genes')
ax2.xaxis.get_label().set_fontsize(40)
ax2.yaxis.get_label().set_fontsize(40)
ax2.grid(False)
ax2.tick_params(labelright=False, right=False)
ax2.set_xlim(0,OverallContrstsSS["size"].max())
ax2.tick_params(labelsize=30)
sns.barplot(data=OverallContrstsSS[OverallContrstsSS['significant'] == 'Down'], x='size', y='Contrast',color="blue",
ci=False, orient='horizontal', dodge=True, ax=ax1)
ax1.invert_xaxis()
ax1.set_title(' '+'Downregulated Genes', fontsize=40)
ax1.set(xlabel='number of genes')
ax1.xaxis.get_label().set_fontsize(40)
ax1.yaxis.get_label().set_fontsize(40)
ax1.tick_params(labelsize=30)
ax1.grid(False)
ax1.set_xlim(OverallContrstsSS["size"].max(),0)
sns.despine()
plt.legend(bbox_to_anchor= (1.6,.5), frameon=False, title="DEGs: {}".format(output_basename.split("_")[1]), fontsize=30,
title_fontsize=30)
fig.savefig("./figures/NumberOfDEGs-{}_Qval-{}_MinlogFC-{}.png".format(output_basename.split("_")[1], QTHRESH,LOGCFTHRESHOLD), bbox_inches='tight')
No handles with labels found to put in legend.